/*	This program produces Appendix Figure C.8 */

***** Set directories 
local dir_clean 	"~/Dropbox/Retirement gaming/clean"
local dir_output 	"~/Dropbox/Retirement gaming/output/dataverse"


use "`dir_clean'/mainsample_medbcw.dta", clear

local listcontrols = " i.t i.ciiu2_1stobs i.ndep_cat_1stobs i.ndep_cat i.year#i.i.ciiu2_1stobs i.year#i.ndep_cat_1stobs "  

 
** EMPLOYEES SMALL FIRMS -- HOURS WORKED (Figure C.8, Panel A) **
preserve
replace time_bcw=-5 if time_bcw==-6
keep if empl==1 &  small_1stobs==1  
*Normalize hours
sum hrsmonth if time_bcw==-1
local meanh=r(mean)
g Hz=hrsmonth/`meanh'
* dummies
tab time_bcw, gen(time_bcw_dums)
drop time_bcw_dums5 
* Regression
reghdfe Hz  time_bcw_dums* , absorb(`listcontrols') vce(cluster i )
* Sample size
global n=e(N)
global N=e(N_clust1)
* Coefficients and sd
for any beta sd: gen X = .
forvalues d = 1(1)13 {
	capture qui replace beta = _b[time_bcw_dums`d'] if time_bcw==`d'-6
	capture qui replace sd   = _se[time_bcw_dums`d'] if time_bcw==`d'-6
}
replace beta=0 if time_bcw==-1
replace sd=0 if time_bcw==-1
* Plot the results
collapse beta* sd*, by(time_bcw)
gen sd_top = beta + 1.96*sd	
gen sd_bot = beta - 1.96*sd
gen sd_top10 = beta + 1.645*sd	
gen sd_bot10 = beta - 1.645*sd
twoway rcap sd_top sd_bot time_bcw, ///
	cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) xlabel(-5(1)7)  xline(-1, lc(green)) ysc(r(-0.03 .09)) ylabel(-0.03(0.03).09,grid)  ///
	|| rcap sd_top10 sd_bot10 time_bcw, cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) ///
	|| scatter beta time_bcw, mcolor(navy) scheme(s1color) lcolor(navy) lpattern(solid) yline(0, lcolor(gs2)) ///
	title("Reported hours per month, relative to time -1") note("N=${n}, Individuals=${N}") ///
	legend(off) ytitle("Estimated Coefficients") xtitle("Years Relative to Predicted Start of BCW")   
graph export "`dir_output'/figureC8a_hours.png",  replace 	
restore

 
*** EMPLOYEES SMALL FIRMS HIGH ATTACHMENT -- HOURS (Figure C.8, Panel B)
preserve
replace time_bcw=-5 if time_bcw==-6
keep if empl==1 & prop_mempl>=.7 &  small_1stobs==1  & tenure1yrs_at49==1
*Normalize hours
sum hrsmonth if time_bcw==-1
local meanh=r(mean)
g Hz=hrsmonth/`meanh'
* dummies
tab time_bcw, gen(time_bcw_dums)
drop time_bcw_dums5 // drop year before ref age, so everything becomes relative to that
* Regression
reghdfe Hz  time_bcw_dums*  , absorb(`listcontrols') vce(cluster i )
* Sample size
global n=e(N)
global N=e(N_clust1)
* Coefficients and sd
for any beta sd: gen X = .
forvalues d = 1(1)13 {
	capture qui replace beta = _b[time_bcw_dums`d'] if time_bcw==`d'-6
	capture qui replace sd   = _se[time_bcw_dums`d'] if time_bcw==`d'-6
}
replace beta=0 if time_bcw==-1
replace sd=0 if time_bcw==-1
* Plot the results
collapse beta* sd*, by(time_bcw)
gen sd_top = beta + 1.96*sd	
gen sd_bot = beta - 1.96*sd
gen sd_top10 = beta + 1.645*sd	
gen sd_bot10 = beta - 1.645*sd
twoway rcap sd_top sd_bot time_bcw, ///
	cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) xlabel(-5(1)7)  xline(-1, lc(green)) ysc(r(-0.03 .09)) ylabel(-0.03(0.03).09,grid) ///
	|| rcap sd_top10 sd_bot10 time_bcw, cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) ///
	|| scatter beta time_bcw, mcolor(navy) scheme(s1color) lcolor(navy) lpattern(solid) yline(0, lcolor(gs2)) ///
	title("Reported hours per month, relative to time -1")  note("N=${n}, Individuals=${N}") ///
	legend(off) ytitle("Estimated Coefficients") xtitle("Years Relative to Predicted Start of BCW")   
graph export "`dir_output'/figureC8b_hours.png",  replace 	
restore



** EMPLOYEES LARGE FIRM -- HOURS (Figure C.8, Panel C) **
preserve
replace time_bcw=-5 if time_bcw==-6
keep if empl==1 & large_1stobs==1
*Normalize wages
sum hrsmonth if time_bcw==-1
local meanh=r(mean)
g Hz=hrsmonth/`meanh'
* dummies
tab time_bcw, gen(time_bcw_dums)
drop time_bcw_dums5 // drop year before ref age, so everything becomes relative to that
* Regression
reghdfe Hz time_bcw_dums*, absorb(`listcontrols') vce(cluster i )
* Sample size
global n=e(N)
global N=e(N_clust1)
* Coefficients and sd
for any beta sd: gen X = .
forvalues d = 1(1)13 {
	capture qui replace beta = _b[time_bcw_dums`d'] if time_bcw==`d'-6
	capture qui replace sd   = _se[time_bcw_dums`d'] if time_bcw==`d'-6
}
replace beta=0 if time_bcw==-1
replace sd=0 if time_bcw==-1
* Plot the results
collapse beta* sd*, by(time_bcw)
gen sd_top = beta + 1.96*sd	
gen sd_bot = beta - 1.96*sd
gen sd_top10 = beta + 1.645*sd	
gen sd_bot10 = beta - 1.645*sd
twoway rcap sd_top sd_bot time_bcw, ///
	cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) xlabel(-5(1)7)  xline(-1, lc(green)) ysc(r(-0.03 .09)) ylabel(-0.03(0.03).09,grid) ///
	|| rcap sd_top10 sd_bot10 time_bcw, cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) ///
	|| scatter beta time_bcw, mcolor(navy) scheme(s1color) lcolor(navy) lpattern(solid) yline(0, lcolor(gs2)) ///
	title("Reported hours per month, relative to time -1") note("N=${n}, Individuals=${N}") ///
	legend(off) ytitle("Estimated Coefficients") xtitle("Years Relative to Predicted Start of BCW")   
graph export "`dir_output'/figureC8c_hours.png",  replace 	
restore


** EMPLOYEES SMALL FIRMS -- wage per hour (Figure C.8, Panel A) **
preserve
replace time_bcw=-5 if time_bcw==-6
keep if empl==1 &  small_1stobs==1  
*Normalize wage per hour
sum wagephr if time_bcw==-1
local meanw=r(mean)
g WPHz=wagephr/`meanw'
* dummies
tab time_bcw, gen(time_bcw_dums)
drop time_bcw_dums5 
* Regression
reghdfe WPHz  time_bcw_dums* , absorb(`listcontrols') vce(cluster i )
* Sample size
global n=e(N)
global N=e(N_clust1)
* Coefficients and sd
for any beta sd: gen X = .
forvalues d = 1(1)13 {
	capture qui replace beta = _b[time_bcw_dums`d'] if time_bcw==`d'-6
	capture qui replace sd   = _se[time_bcw_dums`d'] if time_bcw==`d'-6
}
replace beta=0 if time_bcw==-1
replace sd=0 if time_bcw==-1
* Plot the results
collapse beta* sd*, by(time_bcw)
gen sd_top = beta + 1.96*sd	
gen sd_bot = beta - 1.96*sd
gen sd_top10 = beta + 1.645*sd	
gen sd_bot10 = beta - 1.645*sd
twoway rcap sd_top sd_bot time_bcw, ///
	cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) xlabel(-5(1)7)  xline(-1, lc(green)) ysc(r(-0.03 .09)) ylabel(-0.03(0.03).09,grid)  ///
	|| rcap sd_top10 sd_bot10 time_bcw, cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) ///
	|| scatter beta time_bcw, mcolor(navy) scheme(s1color) lcolor(navy) lpattern(solid) yline(0, lcolor(gs2)) ///
	title("Wage per hour, relative to time -1") note("N=${n}, Individuals=${N}") ///
	legend(off) ytitle("Estimated Coefficients") xtitle("Years Relative to Predicted Start of BCW")   
graph export "`dir_output'/figureC8a_wageperhr.png",  replace 	
restore


*** EMPLOYEES SMALL FIRMS HIGH ATTACHMENT -- wage per hour (Table C.8, Panel B) **
preserve
replace time_bcw=-5 if time_bcw==-6
keep if empl==1 & prop_mempl>=.7 &  small_1stobs==1  & tenure1yrs_at49==1
*Normalize hours
sum wagephr if time_bcw==-1
local meanw=r(mean)
g WPHz=wagephr/`meanw'
* dummies
tab time_bcw, gen(time_bcw_dums)
drop time_bcw_dums5 // drop year before ref age, so everything becomes relative to that
* Regression
reghdfe WPHz  time_bcw_dums*  , absorb(`listcontrols') vce(cluster i )
* Sample size
global n=e(N)
global N=e(N_clust1)
* Coefficients and sd
for any beta sd: gen X = .
forvalues d = 1(1)13 {
	capture qui replace beta = _b[time_bcw_dums`d'] if time_bcw==`d'-6
	capture qui replace sd   = _se[time_bcw_dums`d'] if time_bcw==`d'-6
}
replace beta=0 if time_bcw==-1
replace sd=0 if time_bcw==-1
* Plot the results
collapse beta* sd*, by(time_bcw)
gen sd_top = beta + 1.96*sd	
gen sd_bot = beta - 1.96*sd
gen sd_top10 = beta + 1.645*sd	
gen sd_bot10 = beta - 1.645*sd
twoway rcap sd_top sd_bot time_bcw, ///
	cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) xlabel(-5(1)7)  xline(-1, lc(green)) ysc(r(-0.05 .1)) ylabel(-0.05(0.05).15,grid)  ///
	|| rcap sd_top10 sd_bot10 time_bcw, cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) ///
	|| scatter beta time_bcw, mcolor(navy) scheme(s1color) lcolor(navy) lpattern(solid) yline(0, lcolor(gs2)) ///
	title("Wage per hour, relative to time -1")  note("N=${n}, Individuals=${N}") ///
	legend(off) ytitle("Estimated Coefficients") xtitle("Years Relative to Predicted Start of BCW")   
graph export "`dir_output'/figureC8b_wageperhr.png",  replace 	
restore


*** EMPLOYEES LARGE FIRMS -- wage per hour (Figure C.8, Panel C)
preserve
replace time_bcw=-5 if time_bcw==-6
keep if empl==1 & large_1stobs==1
*Normalize wage per hour
sum wagephr if time_bcw==-1
local meanw=r(mean)
g WPHz=wagephr/`meanw'
* dummies
tab time_bcw, gen(time_bcw_dums)
drop time_bcw_dums5 
* Regression
reghdfe WPHz  time_bcw_dums*, absorb(`listcontrols') vce(cluster i )
* Sample size
global n=e(N)
global N=e(N_clust1)
* Coefficients and sd
for any beta sd: gen X = .
forvalues d = 1(1)13 {
	capture qui replace beta = _b[time_bcw_dums`d'] if time_bcw==`d'-6
	capture qui replace sd   = _se[time_bcw_dums`d'] if time_bcw==`d'-6
}
replace beta=0 if time_bcw==-1
replace sd=0 if time_bcw==-1
* Plot the results
collapse beta* sd*, by(time_bcw)
gen sd_top = beta + 1.96*sd	
gen sd_bot = beta - 1.96*sd
gen sd_top10 = beta + 1.645*sd	
gen sd_bot10 = beta - 1.645*sd
twoway rcap sd_top sd_bot time_bcw, ///
	cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) xlabel(-5(1)7)  xline(-1, lc(green)) ysc(r(-0.03 .09)) ylabel(-0.03(0.03).09,grid) ///
	|| rcap sd_top10 sd_bot10 time_bcw, cmissing(n) lwidth(thin) lcolor(navy)  xsc(r(-5 7)) ///
	|| scatter beta time_bcw, mcolor(navy) scheme(s1color) lcolor(navy) lpattern(solid) yline(0, lcolor(gs2)) ///
	title("Wage per hour, relative to time -1") note("N=${n}, Individuals=${N}") ///
	legend(off) ytitle("Estimated Coefficients") xtitle("Years Relative to Predicted Start of BCW")   
graph export "`dir_output'/figureC8c_wageperhr.png",  replace 	
restore


clear all
exit
 
